Github repo

knitr::opts_chunk$set(message = FALSE, error = FALSE)

library(tidyverse)
library(rvest)
library(magrittr)
library(lubridate)

Understanding our data sources

John Hopkins CSSE

Daily confirmed cases, deaths and recoveries broken down by country/region and state/province, when state/province is available.

confirmed_ts <- read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_19-covid-Confirmed.csv")) %>%
  gather(date, confirmed, -`Province/State`, -`Country/Region`, -Lat, -Long) %>%
  mutate(date = mdy(date))
deaths_ts <- read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_19-covid-Deaths.csv")) %>%
  gather(date, deaths, -`Province/State`, -`Country/Region`, -Lat, -Long) %>%
  mutate(date = mdy(date))
recovered_ts <- read_csv(url("https://raw.githubusercontent.com/CSSEGISandData/COVID-19/master/csse_covid_19_data/csse_covid_19_time_series/time_series_19-covid-Recovered.csv")) %>%
  gather(date, recovered, -`Province/State`, -`Country/Region`, -Lat, -Long) %>%
  mutate(date = mdy(date))

joined <- confirmed_ts %>%
  left_join(deaths_ts) %>%
  left_join(recovered_ts) %>%
  gather(stat, value, confirmed, deaths, recovered) %>%
  arrange(`Country/Region`, `Province/State`, date, stat)

joined %>% mutate_if(is.character, as.factor) %>% summary()
          Province/State         Country/Region       Lat        
 Diamond Princess:  330   US            :40755   Min.   :-41.45  
 Grand Princess  :  330   China         : 5445   1st Qu.: 27.61  
 Adams, IN       :  165   Canada        : 1815   Median : 37.81  
 Alabama         :  165   Australia     : 1485   Mean   : 31.81  
 Alachua, FL     :  165   France        : 1155   3rd Qu.: 42.33  
 (Other)         :50985   United Kingdom:  660   Max.   : 71.71  
 NA's            :24255   (Other)       :25080                   
      Long              date                   stat           value        
 Min.   :-157.86   Min.   :2020-01-22   confirmed:25465   Min.   :    0.0  
 1st Qu.: -93.06   1st Qu.:2020-02-04   deaths   :25465   1st Qu.:    0.0  
 Median : -74.30   Median :2020-02-18   recovered:25465   Median :    0.0  
 Mean   : -36.13   Mean   :2020-02-18                     Mean   :   66.8  
 3rd Qu.:  20.94   3rd Qu.:2020-03-03                     3rd Qu.:    0.0  
 Max.   : 174.89   Max.   :2020-03-16                     Max.   :67798.0  
                                                                           
#' Prepares data.frame for plotting by aggregating, weighting, reordering, 
#' and filtering by location.
#'
#' @param geo_level column name to aggregate (sum) by
#' @param dcr_wts weights to order locations; vector of c(deaths, confirmed, recovered)
#' @param show_top_n show top n locations ordered by dcr_wts
#'
#' @return an aggregated, reordered, and filtered by location data.frame
#'
#' @examples
#' plot_prep("Country/Region", c(100, 1, -1), show_top_n = 25)
#'
#' @export
plot_prep <- function(df, geo_level = "Country/Region", dcr_wts = c(1, 0.01, 0),
                      show_top = 1:10) {
  df %>%
    mutate(location = !!sym(geo_level)) %>%
    # aggregate by location
    group_by(location, date, stat) %>%
    summarise(value = sum(value)) %>%
    ungroup() %>%
    # location order
    group_by(location) %>%
    mutate(
      total_wt =
        value[which(date == max(date) & stat == "deaths")] * dcr_wts[1] +
        value[which(date == max(date) & stat == "confirmed")] * dcr_wts[2] +
        value[which(date == max(date) & stat == "recovered")] * dcr_wts[3]
    ) %>%
    ungroup() %>%
    mutate(total_rank = dense_rank(desc(total_wt))) %>%
    # filter top n
    filter(total_rank %in% show_top) %>%
    # order by rank
    mutate(location = fct_reorder(location, total_rank)) 
}

Cross-country comparison

#' Convenience function to plot country comparisons
plot_country_comps <- function(df) {
  df %>% 
    ggplot(aes(date, value, color=location)) +
    geom_line() + facet_wrap(vars(stat), ncol=1, scales = "free_y") +
    theme(legend.position = "right", legend.title = element_blank()) +
    xlab(element_blank())
}

Top countries

joined %>%
  plot_prep("Country/Region") %>%
  plot_country_comps()

excl China

joined %>%
  filter(!`Country/Region` %in%  c("China")) %>%
  plot_prep("Country/Region") %>%
  plot_country_comps()

excl Italy, Iran, South Korea

joined %>%
  filter(!`Country/Region` %in%  c("China", "Italy", "Iran", "Korea, South")) %>%
  plot_prep("Country/Region") %>%
  plot_country_comps()

excl Spain, France

joined %>%
  filter(!`Country/Region` %in%  c("China", "Italy", "Iran", "Korea, South",
                                   "Spain", "France")) %>%
  plot_prep("Country/Region") %>%
  plot_country_comps()

Comparing growth

#' Convenience function to plot growth comparisons
plot_growth_comps <- function(df) {
  df %>%
    ggplot(aes(date, value)) + geom_line() +
    theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
    facet_wrap(
      vars(location, stat), scales="free_y", ncol = 3,
      labeller=label_wrap_gen(multi_line=FALSE)) +
    theme(strip.background = element_blank(), legend.position = "None") +
    xlab(element_blank())
}

By country

joined %>%
  plot_prep() %>%
  plot_growth_comps()

China

joined %>%
  filter(`Country/Region` == "China") %>%
  plot_prep("Province/State") %>%
  plot_growth_comps()

US

joined %>%
  filter(`Country/Region` == "US") %>%
  plot_prep("Province/State") %>%
  plot_growth_comps()

Canada

joined %>%
  filter(`Country/Region` == "Canada") %>%
  plot_prep("Province/State") %>%
  plot_growth_comps()

Australia

joined %>%
  filter(`Country/Region` == "Australia") %>%
  plot_prep("Province/State") %>%
  plot_growth_comps()

Total deaths since first death

plot_deaths_since_first <- function(df, days_since) {
  df %>%
    group_by(`Country/Region`) %>%
    mutate(total_deaths = value[which(date == max(date) & stat == "deaths")][1]) %>%
    filter(total_deaths > 0) %>%
    plot_prep("Country/Region") %>%
    arrange(location, stat, date) %>%
    group_by(location) %>%
    mutate(first_death_date = date[which(value != 0 & stat == "deaths")][1]) %>%
    mutate(`Days since first death` = date - first_death_date) %>%
    filter(stat == "deaths" & value > 0 &
             `Days since first death` <= days_since) %>%
    ggplot(aes(`Days since first death`, value, color=location)) + geom_line() +
    theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
    theme(strip.background = element_blank(), legend.position = "right",
          legend.title = element_blank()) +
    ylab("Total deaths")
}

10 days since first death

joined %>%
  plot_deaths_since_first(10)

30 days since first death

joined %>%
  plot_deaths_since_first(30)

90 days since first death

joined %>%
  plot_deaths_since_first(90)

Double days

Double days by country/stat

Top 10 countries

joined %>%
  plot_prep(dcr_wts = c(100, 0.1, 0), show_top = 1:10) %>%
  group_by(location, stat) %>%
  arrange(date) %>%
  mutate(
    `Days to double` = if_else(
      value < 10, as.difftime("", units="days"),
      date - date[sapply(value, function(x) { which.min(abs(x / 2 - value)) })])
  ) %>%
  ggplot(aes(date, `Days to double`, color=stat)) +
  geom_point(alpha = 0.5) +
  geom_smooth() +
  #facet_wrap(vars(location), ncol = 3) + 
  facet_grid(vars(location), vars(stat), scales="free") + 
  theme(legend.position = "bottom", legend.title=element_blank()) +
  xlab(element_blank())

excl China

joined %>%
  filter(`Country/Region` != "China") %>%
  plot_prep(dcr_wts = c(100, 0.1, 0), show_top = 1:10) %>%
  group_by(location, stat) %>%
  arrange(date) %>%
  mutate(
    `Days to double` = if_else(
      value < 10, as.difftime("", units="days"),
      date - date[sapply(value, function(x) { which.min(abs(x / 2 - value)) })])
  ) %>%
  ggplot(aes(date, `Days to double`, color=stat)) +
  geom_point(alpha = 0.5) +
  geom_smooth() +
  #facet_wrap(vars(location), ncol = 3) + 
  facet_grid(vars(location), vars(stat), scales="free") + 
  theme(legend.position = "bottom", legend.title=element_blank()) +
  xlab(element_blank())

Days to double deaths since 10th death

h/t: @loeserjohn

Country/Region level

n_deaths <- 10
min_deaths <- n_deaths * 2
top_n_countries <- 1:20
trunc_to_n <- 2
geo_level = "Country/Region"

# TODO(roboton): Refactor into prep and plot functions
joined %>%
  
  # plot prep at country level
  plot_prep(geo_level, show_top = top_n_countries) %>%
  
  # limit to deaths
  filter(stat == "deaths") %>%
  
  # drop countries with less than min_deaths
  group_by(location) %>%
  filter(value[which.max(date)] >= min_deaths) %>%
  ungroup() %>%
  
  # get the first date with more than or equal to n_deaths
  group_by(location) %>%
  mutate(nth_death_date = date[which(value >= n_deaths) - 1][1]) %>%
  mutate(`Days since nth death` = date - nth_death_date) %>%
  
  # Drop days before nth deaths date
  filter(`Days since nth death` > 0) %>%
  
  # Compute Days to Double
  group_by(location) %>%
  
  mutate(
    # set days to double
    double_idx = sapply(value, FUN=function(x) { max(which(value <= x/2)) }),
    `Days to double` = date - date[double_idx]
  ) %>%
  filter(!is.na(`Days to double`)) %>%
  ## plot
  # truncate to second longest time series
  group_by(location) %>%
  mutate(max_days_since_nth_death = max(`Days since nth death`)) %>%
  ungroup() %>%
  mutate(max_days_rank = dense_rank(desc(max_days_since_nth_death))) %>%
  filter(`Days since nth death` <= max(`Days since nth death`[max_days_rank == trunc_to_n])) %>%
  group_by(location) %>%
  mutate(time_points = n()) %>%
  ungroup() %>%
  filter(time_points > 2) %>%
  ggplot(aes(`Days since nth death`, `Days to double`, color=location)) +
  geom_point(alpha = 0.2) +
  geom_smooth(method="auto", se=F) +
  xlab(paste("days since deaths =", n_deaths)) -> p
  ggplotly(p)

State/Province level

n_deaths <- 10
min_deaths <- n_deaths * 2
top_n_countries <- 1:15
trunc_to_n <- 1
geo_level = "Province/State"

# TODO(roboton): Refactor into prep and plot functions
joined %>%
  filter(!is.na(`Province/State`) & `Province/State` != `Country/Region`) %>%
   
  # plot prep at country level
  plot_prep(geo_level, show_top = top_n_countries) %>%
  
  # limit to deaths
  filter(stat == "deaths") %>%
  
  # drop countries with less than min_deaths
  group_by(location) %>%
  filter(value[which.max(date)] >= min_deaths) %>%
  ungroup() %>%
  
  # get the first date with more than or equal to n_deaths
  group_by(location) %>%
  mutate(nth_death_date = date[which(value >= n_deaths) - 1][1]) %>%
  mutate(`Days since nth death` = date - nth_death_date) %>%
  
  # Drop days before nth deaths date
  filter(`Days since nth death` > 0) %>%
  
  # Compute Days to Double
  group_by(location) %>%
  
  mutate(
    # set days to double
    double_idx = sapply(value, FUN=function(x) { max(which(value <= x/2)) }),
    `Days to double` = date - date[double_idx]
  ) %>%
  filter(!is.na(`Days to double`)) %>%
  ## plot
  # truncate to second longest time series
  group_by(location) %>%
  mutate(max_days_since_nth_death = max(`Days since nth death`)) %>%
  ungroup() %>%
  mutate(max_days_rank = dense_rank(desc(max_days_since_nth_death))) %>%
  filter(`Days since nth death` <= max(`Days since nth death`[max_days_rank == trunc_to_n])) %>%
  group_by(location) %>%
  mutate(time_points = n()) %>%
  ungroup() %>%
  filter(time_points > 2) %>%
  ggplot(aes(`Days since nth death`, `Days to double`, color=location)) +
  geom_point(alpha = 0.2) +
  geom_smooth(method="auto", se=F) +
  xlab(paste("days since deaths =", n_deaths)) -> p
  ggplotly(p) 
LS0tCnRpdGxlOiAiQ09WSUQtMTkgRXhwbG9yYXRvcnkiCm91dHB1dDoKICBodG1sX25vdGVib29rOgogICAgY29kZV9mb2xkaW5nOiBoaWRlCiAgICB0aGVtZTogbHVtZW4KICAgIHRvYzogdHJ1ZQogICAgdG9jX2Zsb2F0OiB0cnVlCi0tLQoKW0dpdGh1YiByZXBvXShodHRwczovL2dpdGh1Yi5jb20vcm9ib3Rvbi9jb3ZpZC0xOV9tZXRhKQoKYGBge3Igc2V0dXB9CmtuaXRyOjpvcHRzX2NodW5rJHNldChtZXNzYWdlID0gRkFMU0UsIGVycm9yID0gRkFMU0UpCgpsaWJyYXJ5KHRpZHl2ZXJzZSkKbGlicmFyeShydmVzdCkKbGlicmFyeShtYWdyaXR0cikKbGlicmFyeShsdWJyaWRhdGUpCmBgYAoKIyBVbmRlcnN0YW5kaW5nIG91ciBkYXRhIHNvdXJjZXMKCiMjIEpvaG4gSG9wa2lucyBDU1NFCgpEYWlseSBjb25maXJtZWQgY2FzZXMsIGRlYXRocyBhbmQgcmVjb3ZlcmllcyBicm9rZW4gZG93biBieSBjb3VudHJ5L3JlZ2lvbiBhbmQgc3RhdGUvcHJvdmluY2UsIHdoZW4gc3RhdGUvcHJvdmluY2UgaXMgYXZhaWxhYmxlLgoKYGBge3J9CmNvbmZpcm1lZF90cyA8LSByZWFkX2Nzdih1cmwoImh0dHBzOi8vcmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbS9DU1NFR0lTYW5kRGF0YS9DT1ZJRC0xOS9tYXN0ZXIvY3NzZV9jb3ZpZF8xOV9kYXRhL2Nzc2VfY292aWRfMTlfdGltZV9zZXJpZXMvdGltZV9zZXJpZXNfMTktY292aWQtQ29uZmlybWVkLmNzdiIpKSAlPiUKICBnYXRoZXIoZGF0ZSwgY29uZmlybWVkLCAtYFByb3ZpbmNlL1N0YXRlYCwgLWBDb3VudHJ5L1JlZ2lvbmAsIC1MYXQsIC1Mb25nKSAlPiUKICBtdXRhdGUoZGF0ZSA9IG1keShkYXRlKSkKZGVhdGhzX3RzIDwtIHJlYWRfY3N2KHVybCgiaHR0cHM6Ly9yYXcuZ2l0aHVidXNlcmNvbnRlbnQuY29tL0NTU0VHSVNhbmREYXRhL0NPVklELTE5L21hc3Rlci9jc3NlX2NvdmlkXzE5X2RhdGEvY3NzZV9jb3ZpZF8xOV90aW1lX3Nlcmllcy90aW1lX3Nlcmllc18xOS1jb3ZpZC1EZWF0aHMuY3N2IikpICU+JQogIGdhdGhlcihkYXRlLCBkZWF0aHMsIC1gUHJvdmluY2UvU3RhdGVgLCAtYENvdW50cnkvUmVnaW9uYCwgLUxhdCwgLUxvbmcpICU+JQogIG11dGF0ZShkYXRlID0gbWR5KGRhdGUpKQpyZWNvdmVyZWRfdHMgPC0gcmVhZF9jc3YodXJsKCJodHRwczovL3Jhdy5naXRodWJ1c2VyY29udGVudC5jb20vQ1NTRUdJU2FuZERhdGEvQ09WSUQtMTkvbWFzdGVyL2Nzc2VfY292aWRfMTlfZGF0YS9jc3NlX2NvdmlkXzE5X3RpbWVfc2VyaWVzL3RpbWVfc2VyaWVzXzE5LWNvdmlkLVJlY292ZXJlZC5jc3YiKSkgJT4lCiAgZ2F0aGVyKGRhdGUsIHJlY292ZXJlZCwgLWBQcm92aW5jZS9TdGF0ZWAsIC1gQ291bnRyeS9SZWdpb25gLCAtTGF0LCAtTG9uZykgJT4lCiAgbXV0YXRlKGRhdGUgPSBtZHkoZGF0ZSkpCgpqb2luZWQgPC0gY29uZmlybWVkX3RzICU+JQogIGxlZnRfam9pbihkZWF0aHNfdHMpICU+JQogIGxlZnRfam9pbihyZWNvdmVyZWRfdHMpICU+JQogIGdhdGhlcihzdGF0LCB2YWx1ZSwgY29uZmlybWVkLCBkZWF0aHMsIHJlY292ZXJlZCkgJT4lCiAgYXJyYW5nZShgQ291bnRyeS9SZWdpb25gLCBgUHJvdmluY2UvU3RhdGVgLCBkYXRlLCBzdGF0KQoKam9pbmVkICU+JSBtdXRhdGVfaWYoaXMuY2hhcmFjdGVyLCBhcy5mYWN0b3IpICU+JSBzdW1tYXJ5KCkKYGBgCmBgYHtyfQojJyBQcmVwYXJlcyBkYXRhLmZyYW1lIGZvciBwbG90dGluZyBieSBhZ2dyZWdhdGluZywgd2VpZ2h0aW5nLCByZW9yZGVyaW5nLCAKIycgYW5kIGZpbHRlcmluZyBieSBsb2NhdGlvbi4KIycKIycgQHBhcmFtIGdlb19sZXZlbCBjb2x1bW4gbmFtZSB0byBhZ2dyZWdhdGUgKHN1bSkgYnkKIycgQHBhcmFtIGRjcl93dHMgd2VpZ2h0cyB0byBvcmRlciBsb2NhdGlvbnM7IHZlY3RvciBvZiBjKGRlYXRocywgY29uZmlybWVkLCByZWNvdmVyZWQpCiMnIEBwYXJhbSBzaG93X3RvcF9uIHNob3cgdG9wIG4gbG9jYXRpb25zIG9yZGVyZWQgYnkgZGNyX3d0cwojJwojJyBAcmV0dXJuIGFuIGFnZ3JlZ2F0ZWQsIHJlb3JkZXJlZCwgYW5kIGZpbHRlcmVkIGJ5IGxvY2F0aW9uIGRhdGEuZnJhbWUKIycKIycgQGV4YW1wbGVzCiMnIHBsb3RfcHJlcCgiQ291bnRyeS9SZWdpb24iLCBjKDEwMCwgMSwgLTEpLCBzaG93X3RvcF9uID0gMjUpCiMnCiMnIEBleHBvcnQKcGxvdF9wcmVwIDwtIGZ1bmN0aW9uKGRmLCBnZW9fbGV2ZWwgPSAiQ291bnRyeS9SZWdpb24iLCBkY3Jfd3RzID0gYygxLCAwLjAxLCAwKSwKICAgICAgICAgICAgICAgICAgICAgIHNob3dfdG9wID0gMToxMCkgewogIGRmICU+JQogICAgbXV0YXRlKGxvY2F0aW9uID0gISFzeW0oZ2VvX2xldmVsKSkgJT4lCiAgICAjIGFnZ3JlZ2F0ZSBieSBsb2NhdGlvbgogICAgZ3JvdXBfYnkobG9jYXRpb24sIGRhdGUsIHN0YXQpICU+JQogICAgc3VtbWFyaXNlKHZhbHVlID0gc3VtKHZhbHVlKSkgJT4lCiAgICB1bmdyb3VwKCkgJT4lCiAgICAjIGxvY2F0aW9uIG9yZGVyCiAgICBncm91cF9ieShsb2NhdGlvbikgJT4lCiAgICBtdXRhdGUoCiAgICAgIHRvdGFsX3d0ID0KICAgICAgICB2YWx1ZVt3aGljaChkYXRlID09IG1heChkYXRlKSAmIHN0YXQgPT0gImRlYXRocyIpXSAqIGRjcl93dHNbMV0gKwogICAgICAgIHZhbHVlW3doaWNoKGRhdGUgPT0gbWF4KGRhdGUpICYgc3RhdCA9PSAiY29uZmlybWVkIildICogZGNyX3d0c1syXSArCiAgICAgICAgdmFsdWVbd2hpY2goZGF0ZSA9PSBtYXgoZGF0ZSkgJiBzdGF0ID09ICJyZWNvdmVyZWQiKV0gKiBkY3Jfd3RzWzNdCiAgICApICU+JQogICAgdW5ncm91cCgpICU+JQogICAgbXV0YXRlKHRvdGFsX3JhbmsgPSBkZW5zZV9yYW5rKGRlc2ModG90YWxfd3QpKSkgJT4lCiAgICAjIGZpbHRlciB0b3AgbgogICAgZmlsdGVyKHRvdGFsX3JhbmsgJWluJSBzaG93X3RvcCkgJT4lCiAgICAjIG9yZGVyIGJ5IHJhbmsKICAgIG11dGF0ZShsb2NhdGlvbiA9IGZjdF9yZW9yZGVyKGxvY2F0aW9uLCB0b3RhbF9yYW5rKSkgCn0KYGBgCgojIENyb3NzLWNvdW50cnkgY29tcGFyaXNvbiB7LnRhYnNldCAudGFic2V0LWZhZGUgLnRhYnNldC1waWxsc30KCmBgYHtyfQojJyBDb252ZW5pZW5jZSBmdW5jdGlvbiB0byBwbG90IGNvdW50cnkgY29tcGFyaXNvbnMKcGxvdF9jb3VudHJ5X2NvbXBzIDwtIGZ1bmN0aW9uKGRmKSB7CiAgZGYgJT4lIAogICAgZ2dwbG90KGFlcyhkYXRlLCB2YWx1ZSwgY29sb3I9bG9jYXRpb24pKSArCiAgICBnZW9tX2xpbmUoKSArIGZhY2V0X3dyYXAodmFycyhzdGF0KSwgbmNvbD0xLCBzY2FsZXMgPSAiZnJlZV95IikgKwogICAgdGhlbWUobGVnZW5kLnBvc2l0aW9uID0gInJpZ2h0IiwgbGVnZW5kLnRpdGxlID0gZWxlbWVudF9ibGFuaygpKSArCiAgICB4bGFiKGVsZW1lbnRfYmxhbmsoKSkKfQpgYGAKCgojIyBUb3AgY291bnRyaWVzCgpgYGB7ciBmaWcuaGVpZ2h0ID0gOH0Kam9pbmVkICU+JQogIHBsb3RfcHJlcCgiQ291bnRyeS9SZWdpb24iKSAlPiUKICBwbG90X2NvdW50cnlfY29tcHMoKQpgYGAKCiMjIGV4Y2wgQ2hpbmEKCmBgYHtyIGZpZy5oZWlnaHQgPSA4fQpqb2luZWQgJT4lCiAgZmlsdGVyKCFgQ291bnRyeS9SZWdpb25gICVpbiUgIGMoIkNoaW5hIikpICU+JQogIHBsb3RfcHJlcCgiQ291bnRyeS9SZWdpb24iKSAlPiUKICBwbG90X2NvdW50cnlfY29tcHMoKQpgYGAKCiMjIGV4Y2wgSXRhbHksIElyYW4sIFNvdXRoIEtvcmVhCgpgYGB7ciBmaWcuaGVpZ2h0ID0gOH0Kam9pbmVkICU+JQogIGZpbHRlcighYENvdW50cnkvUmVnaW9uYCAlaW4lICBjKCJDaGluYSIsICJJdGFseSIsICJJcmFuIiwgIktvcmVhLCBTb3V0aCIpKSAlPiUKICBwbG90X3ByZXAoIkNvdW50cnkvUmVnaW9uIikgJT4lCiAgcGxvdF9jb3VudHJ5X2NvbXBzKCkKYGBgCgojIyBleGNsIFNwYWluLCBGcmFuY2UKCmBgYHtyIGZpZy5oZWlnaHQgPSA4fQpqb2luZWQgJT4lCiAgZmlsdGVyKCFgQ291bnRyeS9SZWdpb25gICVpbiUgIGMoIkNoaW5hIiwgIkl0YWx5IiwgIklyYW4iLCAiS29yZWEsIFNvdXRoIiwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAiU3BhaW4iLCAiRnJhbmNlIikpICU+JQogIHBsb3RfcHJlcCgiQ291bnRyeS9SZWdpb24iKSAlPiUKICBwbG90X2NvdW50cnlfY29tcHMoKQpgYGAKCiMgQ29tcGFyaW5nIGdyb3d0aCB7LnRhYnNldCAudGFic2V0LWZhZGUgLnRhYnNldC1waWxsc30KCmBgYHtyfQojJyBDb252ZW5pZW5jZSBmdW5jdGlvbiB0byBwbG90IGdyb3d0aCBjb21wYXJpc29ucwpwbG90X2dyb3d0aF9jb21wcyA8LSBmdW5jdGlvbihkZikgewogIGRmICU+JQogICAgZ2dwbG90KGFlcyhkYXRlLCB2YWx1ZSkpICsgZ2VvbV9saW5lKCkgKwogICAgdGhlbWUoYXhpcy50ZXh0LnggPSBlbGVtZW50X3RleHQoYW5nbGUgPSA0NSwgaGp1c3QgPSAxLCB2anVzdCA9IDEpKSArCiAgICBmYWNldF93cmFwKAogICAgICB2YXJzKGxvY2F0aW9uLCBzdGF0KSwgc2NhbGVzPSJmcmVlX3kiLCBuY29sID0gMywKICAgICAgbGFiZWxsZXI9bGFiZWxfd3JhcF9nZW4obXVsdGlfbGluZT1GQUxTRSkpICsKICAgIHRoZW1lKHN0cmlwLmJhY2tncm91bmQgPSBlbGVtZW50X2JsYW5rKCksIGxlZ2VuZC5wb3NpdGlvbiA9ICJOb25lIikgKwogICAgeGxhYihlbGVtZW50X2JsYW5rKCkpCn0KYGBgCgojIyBCeSBjb3VudHJ5CgpgYGB7ciwgZmlnLndpZHRoPTksIGZpZy5oZWlnaHQ9MTV9CmpvaW5lZCAlPiUKICBwbG90X3ByZXAoKSAlPiUKICBwbG90X2dyb3d0aF9jb21wcygpCmBgYAoKIyMgQ2hpbmEKCmBgYHtyLCBmaWcud2lkdGg9OSwgZmlnLmhlaWdodD0xNX0Kam9pbmVkICU+JQogIGZpbHRlcihgQ291bnRyeS9SZWdpb25gID09ICJDaGluYSIpICU+JQogIHBsb3RfcHJlcCgiUHJvdmluY2UvU3RhdGUiKSAlPiUKICBwbG90X2dyb3d0aF9jb21wcygpCmBgYAoKIyMgVVMKCmBgYHtyLCBmaWcud2lkdGg9OSwgZmlnLmhlaWdodD0xOH0Kam9pbmVkICU+JQogIGZpbHRlcihgQ291bnRyeS9SZWdpb25gID09ICJVUyIpICU+JQogIHBsb3RfcHJlcCgiUHJvdmluY2UvU3RhdGUiKSAlPiUKICBwbG90X2dyb3d0aF9jb21wcygpCmBgYAoKIyMgQ2FuYWRhCgpgYGB7ciwgZmlnLndpZHRoPTksIGZpZy5oZWlnaHQ9MTV9CmpvaW5lZCAlPiUKICBmaWx0ZXIoYENvdW50cnkvUmVnaW9uYCA9PSAiQ2FuYWRhIikgJT4lCiAgcGxvdF9wcmVwKCJQcm92aW5jZS9TdGF0ZSIpICU+JQogIHBsb3RfZ3Jvd3RoX2NvbXBzKCkKYGBgCgojIyBBdXN0cmFsaWEKCmBgYHtyLCBmaWcud2lkdGg9OSwgZmlnLmhlaWdodD0xNX0Kam9pbmVkICU+JQogIGZpbHRlcihgQ291bnRyeS9SZWdpb25gID09ICJBdXN0cmFsaWEiKSAlPiUKICBwbG90X3ByZXAoIlByb3ZpbmNlL1N0YXRlIikgJT4lCiAgcGxvdF9ncm93dGhfY29tcHMoKQpgYGAKCiMgVG90YWwgZGVhdGhzIHNpbmNlIGZpcnN0IGRlYXRoIHsudGFic2V0IC50YWJzZXQtZmFkZSAudGFic2V0LXBpbGxzfQpgYGB7cn0KcGxvdF9kZWF0aHNfc2luY2VfZmlyc3QgPC0gZnVuY3Rpb24oZGYsIGRheXNfc2luY2UpIHsKICBkZiAlPiUKICAgIGdyb3VwX2J5KGBDb3VudHJ5L1JlZ2lvbmApICU+JQogICAgbXV0YXRlKHRvdGFsX2RlYXRocyA9IHZhbHVlW3doaWNoKGRhdGUgPT0gbWF4KGRhdGUpICYgc3RhdCA9PSAiZGVhdGhzIildWzFdKSAlPiUKICAgIGZpbHRlcih0b3RhbF9kZWF0aHMgPiAwKSAlPiUKICAgIHBsb3RfcHJlcCgiQ291bnRyeS9SZWdpb24iKSAlPiUKICAgIGFycmFuZ2UobG9jYXRpb24sIHN0YXQsIGRhdGUpICU+JQogICAgZ3JvdXBfYnkobG9jYXRpb24pICU+JQogICAgbXV0YXRlKGZpcnN0X2RlYXRoX2RhdGUgPSBkYXRlW3doaWNoKHZhbHVlICE9IDAgJiBzdGF0ID09ICJkZWF0aHMiKV1bMV0pICU+JQogICAgbXV0YXRlKGBEYXlzIHNpbmNlIGZpcnN0IGRlYXRoYCA9IGRhdGUgLSBmaXJzdF9kZWF0aF9kYXRlKSAlPiUKICAgIGZpbHRlcihzdGF0ID09ICJkZWF0aHMiICYgdmFsdWUgPiAwICYKICAgICAgICAgICAgIGBEYXlzIHNpbmNlIGZpcnN0IGRlYXRoYCA8PSBkYXlzX3NpbmNlKSAlPiUKICAgIGdncGxvdChhZXMoYERheXMgc2luY2UgZmlyc3QgZGVhdGhgLCB2YWx1ZSwgY29sb3I9bG9jYXRpb24pKSArIGdlb21fbGluZSgpICsKICAgIHRoZW1lKGF4aXMudGV4dC54ID0gZWxlbWVudF90ZXh0KGFuZ2xlID0gNDUsIGhqdXN0ID0gMSwgdmp1c3QgPSAxKSkgKwogICAgdGhlbWUoc3RyaXAuYmFja2dyb3VuZCA9IGVsZW1lbnRfYmxhbmsoKSwgbGVnZW5kLnBvc2l0aW9uID0gInJpZ2h0IiwKICAgICAgICAgIGxlZ2VuZC50aXRsZSA9IGVsZW1lbnRfYmxhbmsoKSkgKwogICAgeWxhYigiVG90YWwgZGVhdGhzIikKfQpgYGAKCiMjIDEwIGRheXMgc2luY2UgZmlyc3QgZGVhdGgKYGBge3IgZmlnLndpZHRoID0gOH0Kam9pbmVkICU+JQogIHBsb3RfZGVhdGhzX3NpbmNlX2ZpcnN0KDEwKQpgYGAKCiMjIDMwIGRheXMgc2luY2UgZmlyc3QgZGVhdGgKYGBge3IgZmlnLndpZHRoID0gOH0Kam9pbmVkICU+JQogIHBsb3RfZGVhdGhzX3NpbmNlX2ZpcnN0KDMwKQpgYGAKCiMjIDkwIGRheXMgc2luY2UgZmlyc3QgZGVhdGgKYGBge3IgZmlnLndpZHRoID0gOH0Kam9pbmVkICU+JQogIHBsb3RfZGVhdGhzX3NpbmNlX2ZpcnN0KDkwKQpgYGAKCiMgRG91YmxlIGRheXMKCiMjIERvdWJsZSBkYXlzIGJ5IGNvdW50cnkvc3RhdCB7LnRhYnNldCAudGFic2V0LWZhZGUgLnRhYnNldC1waWxsc30KCiMjIyBUb3AgMTAgY291bnRyaWVzCgpgYGB7ciBmaWcuaGVpZ2h0PTEwLCBtZXNzYWdlPUZBTFNFLCB3YXJuaW5nPUZBTFNFfQpqb2luZWQgJT4lCiAgcGxvdF9wcmVwKGRjcl93dHMgPSBjKDEwMCwgMC4xLCAwKSwgc2hvd190b3AgPSAxOjEwKSAlPiUKICBncm91cF9ieShsb2NhdGlvbiwgc3RhdCkgJT4lCiAgYXJyYW5nZShkYXRlKSAlPiUKICBtdXRhdGUoCiAgICBgRGF5cyB0byBkb3VibGVgID0gaWZfZWxzZSgKICAgICAgdmFsdWUgPCAxMCwgYXMuZGlmZnRpbWUoIiIsIHVuaXRzPSJkYXlzIiksCiAgICAgIGRhdGUgLSBkYXRlW3NhcHBseSh2YWx1ZSwgZnVuY3Rpb24oeCkgeyB3aGljaC5taW4oYWJzKHggLyAyIC0gdmFsdWUpKSB9KV0pCiAgKSAlPiUKICBnZ3Bsb3QoYWVzKGRhdGUsIGBEYXlzIHRvIGRvdWJsZWAsIGNvbG9yPXN0YXQpKSArCiAgZ2VvbV9wb2ludChhbHBoYSA9IDAuNSkgKwogIGdlb21fc21vb3RoKCkgKwogICNmYWNldF93cmFwKHZhcnMobG9jYXRpb24pLCBuY29sID0gMykgKyAKICBmYWNldF9ncmlkKHZhcnMobG9jYXRpb24pLCB2YXJzKHN0YXQpLCBzY2FsZXM9ImZyZWUiKSArIAogIHRoZW1lKGxlZ2VuZC5wb3NpdGlvbiA9ICJib3R0b20iLCBsZWdlbmQudGl0bGU9ZWxlbWVudF9ibGFuaygpKSArCiAgeGxhYihlbGVtZW50X2JsYW5rKCkpCmBgYAoKIyMjIGV4Y2wgQ2hpbmEKCmBgYHtyIGZpZy5oZWlnaHQ9MTAsIG1lc3NhZ2U9RkFMU0UsIHdhcm5pbmc9RkFMU0V9CmpvaW5lZCAlPiUKICBmaWx0ZXIoYENvdW50cnkvUmVnaW9uYCAhPSAiQ2hpbmEiKSAlPiUKICBwbG90X3ByZXAoZGNyX3d0cyA9IGMoMTAwLCAwLjEsIDApLCBzaG93X3RvcCA9IDE6MTApICU+JQogIGdyb3VwX2J5KGxvY2F0aW9uLCBzdGF0KSAlPiUKICBhcnJhbmdlKGRhdGUpICU+JQogIG11dGF0ZSgKICAgIGBEYXlzIHRvIGRvdWJsZWAgPSBpZl9lbHNlKAogICAgICB2YWx1ZSA8IDEwLCBhcy5kaWZmdGltZSgiIiwgdW5pdHM9ImRheXMiKSwKICAgICAgZGF0ZSAtIGRhdGVbc2FwcGx5KHZhbHVlLCBmdW5jdGlvbih4KSB7IHdoaWNoLm1pbihhYnMoeCAvIDIgLSB2YWx1ZSkpIH0pXSkKICApICU+JQogIGdncGxvdChhZXMoZGF0ZSwgYERheXMgdG8gZG91YmxlYCwgY29sb3I9c3RhdCkpICsKICBnZW9tX3BvaW50KGFscGhhID0gMC41KSArCiAgZ2VvbV9zbW9vdGgoKSArCiAgI2ZhY2V0X3dyYXAodmFycyhsb2NhdGlvbiksIG5jb2wgPSAzKSArIAogIGZhY2V0X2dyaWQodmFycyhsb2NhdGlvbiksIHZhcnMoc3RhdCksIHNjYWxlcz0iZnJlZSIpICsgCiAgdGhlbWUobGVnZW5kLnBvc2l0aW9uID0gImJvdHRvbSIsIGxlZ2VuZC50aXRsZT1lbGVtZW50X2JsYW5rKCkpICsKICB4bGFiKGVsZW1lbnRfYmxhbmsoKSkKYGBgCgojIyAgRGF5cyB0byBkb3VibGUgZGVhdGhzIHNpbmNlIDEwdGggZGVhdGggey50YWJzZXQgLnRhYnNldC1mYWRlIC50YWJzZXQtcGlsbHN9CgpoL3Q6IFtcQGxvZXNlcmpvaG5dKGh0dHBzOi8vdHdpdHRlci5jb20vbG9lc2Vyam9obikKCiMjIyBDb3VudHJ5L1JlZ2lvbiBsZXZlbAoKYGBge3IgZmlnLndpZHRoID0gNywgd2FybmluZyA9IEZBTFNFLCBlcnJvciA9IEZBTFNFfQpuX2RlYXRocyA8LSAxMAptaW5fZGVhdGhzIDwtIG5fZGVhdGhzICogMgp0b3Bfbl9jb3VudHJpZXMgPC0gMToyMAp0cnVuY190b19uIDwtIDIKZ2VvX2xldmVsID0gIkNvdW50cnkvUmVnaW9uIgoKIyBUT0RPKHJvYm90b24pOiBSZWZhY3RvciBpbnRvIHByZXAgYW5kIHBsb3QgZnVuY3Rpb25zCmpvaW5lZCAlPiUKICAKICAjIHBsb3QgcHJlcCBhdCBjb3VudHJ5IGxldmVsCiAgcGxvdF9wcmVwKGdlb19sZXZlbCwgc2hvd190b3AgPSB0b3Bfbl9jb3VudHJpZXMpICU+JQogIAogICMgbGltaXQgdG8gZGVhdGhzCiAgZmlsdGVyKHN0YXQgPT0gImRlYXRocyIpICU+JQogIAogICMgZHJvcCBjb3VudHJpZXMgd2l0aCBsZXNzIHRoYW4gbWluX2RlYXRocwogIGdyb3VwX2J5KGxvY2F0aW9uKSAlPiUKICBmaWx0ZXIodmFsdWVbd2hpY2gubWF4KGRhdGUpXSA+PSBtaW5fZGVhdGhzKSAlPiUKICB1bmdyb3VwKCkgJT4lCiAgCiAgIyBnZXQgdGhlIGZpcnN0IGRhdGUgd2l0aCBtb3JlIHRoYW4gb3IgZXF1YWwgdG8gbl9kZWF0aHMKICBncm91cF9ieShsb2NhdGlvbikgJT4lCiAgbXV0YXRlKG50aF9kZWF0aF9kYXRlID0gZGF0ZVt3aGljaCh2YWx1ZSA+PSBuX2RlYXRocykgLSAxXVsxXSkgJT4lCiAgbXV0YXRlKGBEYXlzIHNpbmNlIG50aCBkZWF0aGAgPSBkYXRlIC0gbnRoX2RlYXRoX2RhdGUpICU+JQogIAogICMgRHJvcCBkYXlzIGJlZm9yZSBudGggZGVhdGhzIGRhdGUKICBmaWx0ZXIoYERheXMgc2luY2UgbnRoIGRlYXRoYCA+IDApICU+JQogIAogICMgQ29tcHV0ZSBEYXlzIHRvIERvdWJsZQogIGdyb3VwX2J5KGxvY2F0aW9uKSAlPiUKICAKICBtdXRhdGUoCiAgICAjIHNldCBkYXlzIHRvIGRvdWJsZQogICAgZG91YmxlX2lkeCA9IHNhcHBseSh2YWx1ZSwgRlVOPWZ1bmN0aW9uKHgpIHsgbWF4KHdoaWNoKHZhbHVlIDw9IHgvMikpIH0pLAogICAgYERheXMgdG8gZG91YmxlYCA9IGRhdGUgLSBkYXRlW2RvdWJsZV9pZHhdCiAgKSAlPiUKICBmaWx0ZXIoIWlzLm5hKGBEYXlzIHRvIGRvdWJsZWApKSAlPiUKICAjIyBwbG90CiAgIyB0cnVuY2F0ZSB0byBzZWNvbmQgbG9uZ2VzdCB0aW1lIHNlcmllcwogIGdyb3VwX2J5KGxvY2F0aW9uKSAlPiUKICBtdXRhdGUobWF4X2RheXNfc2luY2VfbnRoX2RlYXRoID0gbWF4KGBEYXlzIHNpbmNlIG50aCBkZWF0aGApKSAlPiUKICB1bmdyb3VwKCkgJT4lCiAgbXV0YXRlKG1heF9kYXlzX3JhbmsgPSBkZW5zZV9yYW5rKGRlc2MobWF4X2RheXNfc2luY2VfbnRoX2RlYXRoKSkpICU+JQogIGZpbHRlcihgRGF5cyBzaW5jZSBudGggZGVhdGhgIDw9IG1heChgRGF5cyBzaW5jZSBudGggZGVhdGhgW21heF9kYXlzX3JhbmsgPT0gdHJ1bmNfdG9fbl0pKSAlPiUKICBncm91cF9ieShsb2NhdGlvbikgJT4lCiAgbXV0YXRlKHRpbWVfcG9pbnRzID0gbigpKSAlPiUKICB1bmdyb3VwKCkgJT4lCiAgZmlsdGVyKHRpbWVfcG9pbnRzID4gMikgJT4lCiAgZ2dwbG90KGFlcyhgRGF5cyBzaW5jZSBudGggZGVhdGhgLCBgRGF5cyB0byBkb3VibGVgLCBjb2xvcj1sb2NhdGlvbikpICsKICBnZW9tX3BvaW50KGFscGhhID0gMC4yKSArCiAgZ2VvbV9zbW9vdGgobWV0aG9kPSJhdXRvIiwgc2U9RikgKwogIHhsYWIocGFzdGUoImRheXMgc2luY2UgZGVhdGhzID0iLCBuX2RlYXRocykpIC0+IHAKICBnZ3Bsb3RseShwKQpgYGAKCiMjIyBTdGF0ZS9Qcm92aW5jZSBsZXZlbAoKYGBge3IgZmlnLndpZHRoID0gNywgd2FybmluZyA9IEZBTFNFLCBlcnJvciA9IEZBTFNFfQpuX2RlYXRocyA8LSAxMAptaW5fZGVhdGhzIDwtIG5fZGVhdGhzICogMgp0b3Bfbl9jb3VudHJpZXMgPC0gMToxNQp0cnVuY190b19uIDwtIDEKZ2VvX2xldmVsID0gIlByb3ZpbmNlL1N0YXRlIgoKIyBUT0RPKHJvYm90b24pOiBSZWZhY3RvciBpbnRvIHByZXAgYW5kIHBsb3QgZnVuY3Rpb25zCmpvaW5lZCAlPiUKICBmaWx0ZXIoIWlzLm5hKGBQcm92aW5jZS9TdGF0ZWApICYgYFByb3ZpbmNlL1N0YXRlYCAhPSBgQ291bnRyeS9SZWdpb25gKSAlPiUKICAgCiAgIyBwbG90IHByZXAgYXQgY291bnRyeSBsZXZlbAogIHBsb3RfcHJlcChnZW9fbGV2ZWwsIHNob3dfdG9wID0gdG9wX25fY291bnRyaWVzKSAlPiUKICAKICAjIGxpbWl0IHRvIGRlYXRocwogIGZpbHRlcihzdGF0ID09ICJkZWF0aHMiKSAlPiUKICAKICAjIGRyb3AgY291bnRyaWVzIHdpdGggbGVzcyB0aGFuIG1pbl9kZWF0aHMKICBncm91cF9ieShsb2NhdGlvbikgJT4lCiAgZmlsdGVyKHZhbHVlW3doaWNoLm1heChkYXRlKV0gPj0gbWluX2RlYXRocykgJT4lCiAgdW5ncm91cCgpICU+JQogIAogICMgZ2V0IHRoZSBmaXJzdCBkYXRlIHdpdGggbW9yZSB0aGFuIG9yIGVxdWFsIHRvIG5fZGVhdGhzCiAgZ3JvdXBfYnkobG9jYXRpb24pICU+JQogIG11dGF0ZShudGhfZGVhdGhfZGF0ZSA9IGRhdGVbd2hpY2godmFsdWUgPj0gbl9kZWF0aHMpIC0gMV1bMV0pICU+JQogIG11dGF0ZShgRGF5cyBzaW5jZSBudGggZGVhdGhgID0gZGF0ZSAtIG50aF9kZWF0aF9kYXRlKSAlPiUKICAKICAjIERyb3AgZGF5cyBiZWZvcmUgbnRoIGRlYXRocyBkYXRlCiAgZmlsdGVyKGBEYXlzIHNpbmNlIG50aCBkZWF0aGAgPiAwKSAlPiUKICAKICAjIENvbXB1dGUgRGF5cyB0byBEb3VibGUKICBncm91cF9ieShsb2NhdGlvbikgJT4lCiAgCiAgbXV0YXRlKAogICAgIyBzZXQgZGF5cyB0byBkb3VibGUKICAgIGRvdWJsZV9pZHggPSBzYXBwbHkodmFsdWUsIEZVTj1mdW5jdGlvbih4KSB7IG1heCh3aGljaCh2YWx1ZSA8PSB4LzIpKSB9KSwKICAgIGBEYXlzIHRvIGRvdWJsZWAgPSBkYXRlIC0gZGF0ZVtkb3VibGVfaWR4XQogICkgJT4lCiAgZmlsdGVyKCFpcy5uYShgRGF5cyB0byBkb3VibGVgKSkgJT4lCiAgIyMgcGxvdAogICMgdHJ1bmNhdGUgdG8gc2Vjb25kIGxvbmdlc3QgdGltZSBzZXJpZXMKICBncm91cF9ieShsb2NhdGlvbikgJT4lCiAgbXV0YXRlKG1heF9kYXlzX3NpbmNlX250aF9kZWF0aCA9IG1heChgRGF5cyBzaW5jZSBudGggZGVhdGhgKSkgJT4lCiAgdW5ncm91cCgpICU+JQogIG11dGF0ZShtYXhfZGF5c19yYW5rID0gZGVuc2VfcmFuayhkZXNjKG1heF9kYXlzX3NpbmNlX250aF9kZWF0aCkpKSAlPiUKICBmaWx0ZXIoYERheXMgc2luY2UgbnRoIGRlYXRoYCA8PSBtYXgoYERheXMgc2luY2UgbnRoIGRlYXRoYFttYXhfZGF5c19yYW5rID09IHRydW5jX3RvX25dKSkgJT4lCiAgZ3JvdXBfYnkobG9jYXRpb24pICU+JQogIG11dGF0ZSh0aW1lX3BvaW50cyA9IG4oKSkgJT4lCiAgdW5ncm91cCgpICU+JQogIGZpbHRlcih0aW1lX3BvaW50cyA+IDIpICU+JQogIGdncGxvdChhZXMoYERheXMgc2luY2UgbnRoIGRlYXRoYCwgYERheXMgdG8gZG91YmxlYCwgY29sb3I9bG9jYXRpb24pKSArCiAgZ2VvbV9wb2ludChhbHBoYSA9IDAuMikgKwogIGdlb21fc21vb3RoKG1ldGhvZD0iYXV0byIsIHNlPUYpICsKICB4bGFiKHBhc3RlKCJkYXlzIHNpbmNlIGRlYXRocyA9Iiwgbl9kZWF0aHMpKSAtPiBwCiAgZ2dwbG90bHkocCkgCmBgYA==